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Abstract.  A model  of  a sunspot  penumbra,  including  (he  elfects  of  magnetic  field,  compressibility,  and 
buoyancy,  is  studied  in  order  to  identify  the  mode  of  running  penumbral  waves.  It  is  found  that  the 
penumbral  waves  may  be  identified  with  gravity-modified  magneto-acoustic  waves  of  the  'plus'  type 
that  are  vertically  trapped  at  photospheric  levels.  Although  most  of  the  wave  energy  is  contained  in 
the  penumbral  photosphere  and  subphotosphere,  the  maximum  vertical  velocity  occurs  in  the  chro- 
mosphere where  ft)  the  waves  are  evanescent  and  (ii)  the  vertical  velocity  is  in  fact  observed  (in  H v). 

1.  Introduction 

Recent  observations  have  disclosed  an  interesting  pattern  of  velocity  fields  in  sun- 
spots. The  most  recent  discovery  is  that  of  waves  propagating  radially  outward  in 
sunspot  penumbrae  (Zirin  and  Stein,  1972;  Giovanelli,  1972).  Zirin  and  Stein  refer  to 
these  waves  as  running  penumbral  waves.  With  the  further  observations  of  Giovanelli 
(1974),  we  now  have  a fairly  clear  picture  of  the  properties  of  these  waves.  The  purpose 
of  this  theoretical  paper  is  to  study  possible  wave  modes  in  a model  of  a sunspot 
penumbra  in  order  to  identify  the  mode  of  the  running  penumbral  waves.  We  shall 
argue  that  the  running  penumbral  waves  are  gravity-modified  magneto-acoustic 
waves  (of  the  ‘plus’  type)  that  are  vertically  trapped  at  photospheric  levels. 

Giovanelli  (1974)  has  summarized  the  observations  of  running  penumbral  waves, 
and  he  presents  the  following  picture.  The  waves  are  observed  in  Hot  by  means  of 
their  line-of-sight  velocity.  They  occur  in  almost  every  sizable  spot  with  a regular 
stable  structure,  but  only  rarely  in  active  spots  with  complex  structure.  The  waves 
travel  outward  in  the  penumbra  at  a typical  speed  of  15  km  s'1.  The  observed  waves 
have  periods  in  the  range  180-240  s and  horizontal  wavelengths  in  the  range  2350- 
3800  km.  Observations  near  the  limb  have  failed  to  reveal  any  horizontal  motions 
associated  with  the  penumbral  waves,  so  the  wave  motion  is  predominantly  vertical 
in  H*. 

Thus  far  no  detailed  theoretical  study  of  the  mode  of  the  running  penumbral  waves 
has  appeared,  although  Moore  (1973)  has  studied  the  related  problem  of  the  genera- 
tion of  penumbral  waves  in  the  umbra.  Zirin  and  Stein  (1972)  tentatively  identified 
the  penumbral  waves  as  sound  waves,  whereas  Giovanelli  (1972,  1974)  identified  them 
as  Alfven  waves.  The  penumbral  waves,  with  their  predominant  vertical  motions,  no 
doubt  involve  the  combined  effects  of  restoring  forces  due  to  compressibility,  magnetic 

* Also  National  Science  Foundation  Predoctoral  Trainee,  Dept,  of  Mechanical  and  Aerospace 
Sciences,  University  of  Rochester,  Rochester.  N.Y..  U.S.A. 
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field,  and  buoyancy,  and  a complete  theory  should  account  for  all  three  effects.  This 
is  done  in  the  present  paper. 

In  studying  penumbral  waves,  we  face  a difficulty,  in  that  there  seems  to  be  no 
complete,  generally  accepted  penumbral  model  on  which  to  base  our  calculations. 
We  have  therefore  constructed  a penumbral  model  for  use  in  studying  wave  modes. 
This  model  (presented  in  Section  3),  while  simple  enough  to  permit  analysis  of  wave 
modes,  nevertheless  reproduces  all  of  the  relevant  features  of  penumbral  structure, 
and  is  in  reasonable  quantitative  agreement  with  observations.  We  have  assumed  the 
penumbral  magnetic  held  to  be  purely  horizontal,  but  varying  with  height.  True  pe- 
numbral magnetic  fields  are  not  purely  horizontal,  although  they  may  be  very  nearly 
so  (Nishi  and  Makita,  1973).  There  is  some  disagreement  over  the  inclination  of  the 
magnetic  field  in  a penumbra  (see  Beckers  and  Schrdter  (1969)  for  a summary  of 
observations).  The  assumption  of  a horizontal  field  here  is  mostly  a matter  of  con- 
venience; the  basic  mechanism  we  propose  for  the  vertical  trapping  of  penumbral 
waves  will  also  work  for  an  inclined  field.  We  have  also  taken  our  model  to  be  horizon- 
tally uniform  - that  is,  we  have  not  tried  to  represent  the  horizontal  filamentary 
structure  of  a penumbra  or  the  radial  geometry. 

In  Section  2 we  present  the  basic  equations  for  waves  in  our  model  penumbral 
atmosphere.  The  basic  atmosphere  is  completely  characterized  in  these  equations  by 
the  vertical  distribution  of  three  parameters:  the  sound  speed  c,  the  Alfven  velocity 
tA,  and  the  local  density  scale  height  H.  In  order  to  illustrate  the  properties  of  the 
various  wave  modes  that  can  occur,  we  study  the  dispersion  relation  that  holds  in  the 
case  of  constant  c,  rA,  and  //.  In  Section  3 we  present  the  basic  penumbral  model  in 
terms  of  the  distributions  of  c,  eA,  and  H with  height.  In  Section  4 we  show  that  the 
penumbral  waves  may  be  identified  with  ‘plus'  modes  that  are  trapped  in  the  photo- 
spheric-subphotospheric  region  in  our  model.  We  discuss  these  modes  further  in 
Section  5. 

2.  Basic  Equations  and  Dispersion  Relations 

In  our  simplified  treatment  of  a sunspot  penumbra  we  shall  ignore  the  radial  spreading 
of  magnetic  field  lines,  and  consider  the  undisturbed  magnetic  field  to  be  purely 
horizontal  (in  the  .v-direction)  and  varying  with  height  z;  i.e.,  B(l  = (B0(z).  0,0). 
We  assume  the  field  permeates  an  inviscid,  perfectly  conducting,  plane  stratified 
atmosphere  with  constant  acceleration  of  gravity  g (=0.274  km  s'2)  in  the  negative 
z-direction.  The  undisturbed  pressure,  density,  and  temperature  are  denoted  by 
/>o(z),  £>o(z).  and  7"0(z),  respectively.  The  atmosphere  is  in  hydrostatic  equilibrium, 
so  that 

d / B20\ 

dzlPo  + 8j  = -^-  (,) 

We  then  consider  small  adiabatic  perturbations  of  this  equilibrium  atmosphere. 
We  consider  wave  vectors  only  in  the  .vz  plane,  and  assume  that  the  perturbation 
velocity  u=(u,  v,  »•)  has  the  form  u = u exp/ (kxx— u>l ),  with  u = u(z)=  (i?(z),  f(z). 
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it (z ) ).  Starting  with  the  linearized  equations  of  continuity,  energy,  and  momentum, 
we  can  eliminate  the  pressure  and  density  perturbations,  and  arrive  at  the  set  of  three 
linearized  momentum  equations  for  the  velocity  components  u.  f,  and  if,  in  which  the 
basic  atmosphere  is  completely  described  by  the  sound  speed  r(z),  the  Alfven  velocity 
rA(z),  and  the  local  density  scale  height  //  (z).  defined  by 


_ / cPo\ 

\<  L'o/s 


Bo  1 = _ 1 de0 

4ng0 ' H c'o 


The  linearized  momentum  equations  are  the  following: 


(to2  — c2k 2)  u + ikxc2  — ikxgw  = 0, 


(w2-vik2x)v  = 0. 


(c2  + V 


,,  d2if  rd  , , _ C'2  + uA”|  dif  , , , 

J dz2  + [dz  (C  + "a)  ~~H~  J dz  + (l°  ~ Va><x)  “ 

f dc2  c2  "1  , dw 


+ i*‘[dz  ~ H +aJ"  + iA"c  dz  = 0'  <S> 

We  can  eliminate  the  horizontal  velocity  components  u and  0 from  the  system  of 
Equations  (3H5)  to  obtain  a single  equation  for  the  vertical  velocity  if',  in  the  form 

d2if'  dif 

, + -d  (z)  — + B(z)  if  = 0,  (6) 

dz  dz 


A (z)  = - + 

’ H 


1 r oy  df2  , , , , dr n 

■ H + iz  + <“  - ' **)  dTj * 

x [eA(to2  - c2k2x)  + cW]  1 


= [V  - M 


.ikl)(a>2-c2kl)-glg-  \ki- 


co2A2  dc2 
(to2  - c2k2x)  dz 


d;2l[ri(W2-f2Ai)  + f2to2]-'. 

dz  _ 


With  the  transformation 


4>(z)  = H'(z)exp 


Equation  (6)  assumes  the  form 


2 + J (kXi  (o;  z)<f>  = 0, 


where 


f=B-\A2  - J 


d4 

dz 


(M) 


Here,  |<^)|2  is  roughly  proportional  to  the  kinetic  energy  of  wave  motion.  Equation 
(10)  is  the  propagation  equation  for  waves  in  an  atmosphere  with  vertically  varying 
c2,  i'a,  and  //.  For  given  distributions  of  c2  (;),  r*  (z),  and  //  (z),  we  can  use  the  expres- 
sion for  / (c)  (Equation  (ID)  to  distinguish  roughly  between  local  regions  where  a 
wave  with  a particular  frequency  10  and  horizontal  wavelength  kx  is  vertically  propa- 
gating (/>0)  or  vertically  evanescent  (/<0).  We  shall  use  this  approach  for  our 
penumbral  model  in  Section  4. 

The  simplest  case  to  study  is  that  of  constant  c2,  (•*,  and  H.  We  shall  consider  this 
case  now  in  order  to  show  the  kinds  of  wave  modes  which  can  occur.  Although  this 
case  does  not  apply  strictly  to  a real  penumbra,  or  in  fact  to  our  penumbral  model,  we 
can  nevertheless  apply  the  resulting  dispersion  relation  locally  to  get  an  approximate 
picture  of  the  wave  modes.  We  can  also  approximate  a continuous  vertical  variation 
of  the  parameters  in  the  penumbra  by  a series  of  layers  in  which  they  are  constant.  In 
the  case  of  constant  c2,  v\,  and  //.  Equations  (3)— (5)  have  constant  coefficients,  and 
we  can  assume  a solution  of  the  form  u(z)  = u e\p(ik.z+z/2H ),  where  u=  (//,  f.  w)  is 
a constant  vector.  Here,  the  factor  exp(z/2//)  accounts  for  the  fact  that,  to  conserve 
energy,  the  perturbation  amplitude  must  grow  as  the  density  decreases.  The  system  of 
Equations  (3)— (5)  then  becomes  (cf.,  Yu,  1965): 


The  dispersion  relation  for  waves  is  obtained  from  the  condition  for  nonzero  solu- 
tions of  the  homogeneous  system  ( 1 2)— ( 1 4),  i.e.,  the  vanishing  of  the  determinant  of 
the  coefficients.  This  yields  the  dispersion  relation 


(w2  - 1>^;)  •!«  - (‘  + i’a)  I k*  + a,  + 4I{1 


id2  + 


-f-  c2v^k2  ( k2  T k2  T , 

\ x\  x . 4H2 


~ 9 (d  ~ ^ j = 0 • (15) 


We  now  discuss  the  various  wave  modes  given  by  this  dispersion  relation,  with  some 
comments  on  their  relation  to  running  penumbral  waves. 
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2.1.  The:  pure  acoustic  mode 

The  dispersion  relation  (15)  is  satisfied  with  w2—c2k\  for  the  particular  imaginary 
value  of  A.  given  by  /A.  = (1/2//)—  (g/c2).  From  Equations  (I2)-(14)  we  see  the  cor- 
responding motion  has  i?#0,  whereas  f = vr  = 0.  This  mode  is  thus  a purely  compres- 
sional,  acoustic  mode  with  dispersion  relation  identical  to  that  of  a homogeneous 
isothermal  gas.  This  can  occur  since  the  motion  is  horizontal  and  parallel  to  the 
magnetic  field,  and  thus  there  is  no  contribution  from  buoyant  or  magnetic  forces. 
The  amplitude  behaves  as 

and  thus  grows  exponentially  with  height  (since  H <c2  Ig  for  stability  of  the  unpertur- 
bed atmosphere  (Yu,  1 965)).  The  total  momentum  and  energy  are  finite,  however, 
provided  there  is  a lower  bound,  ry  confining  the  motion  to  a semi-infinite  range  of  r. 
This  mode  is  identical  to  the  Lamb  mode  in  the  non-magnetic  case,  except  that  here 
the  scale  height  H is  modified  by  the  magnetic  field.  Since  this  mode  has  no  vertical 
motion,  it  can  net  be  associated  with  the  running  penumbral  waxes. 


2.2.  The  pure  alfven  mode 

A root  of  the  dispersion  relation  (15).  for  arbitrary  A.,  is  given  by  oj2  = c2k2x,  which  is 
the  same  as  the  dispersion  relation  for  a pure  Alfven  wave  in  a homogeneous  atmos- 
phere with  uniform  magnetic  field  With  ar  = r*Aj,  however.  the  other  factor  in  the 
dispersion  relation  ( 15)  is  in  genera!  not  zero,  and  thus  Equations  ( 1 2)— ( 1 4)  sho\xj  that 
the  motion  has  f#0,  whereas  i7=m=0.  The  motion  is  purely  horizontal  and  purely 
transverse  to  the  magnetic  field.  Thus,  the  pure  Alfven  mode  in  a stratified  atmosphere 
is.  so  to  speak,  plane  polarized,  with  no  vertical  motions.  With  n=0  there  is  no  con- 
tribution from  the  buoyancy  force,  and.  further,  with  i7  = 0 the  motion  is  incompres- 
sible. so  the  wave  behaves  as  a pure  Alfven  wave  (with  the  amplitude  factor 
exp  (c/2/-/)).  Again,  since  this  mode  has  no  vertical  motions,  it  can  not  be  associated 
with  the  running  penumbral  waves.  We  turn  now  to  the  remaining  roots  of  the  dis- 
persion relation  (15),  which  do  permit  vertical  motions. 


2.3.  The  plus  and  minus  modes 

The  remaining  roots  of  the  dispersion  relation  (15)  are  given  by 


(t) 


2 

± 


= H(c2  + e, 


2'  ‘ kl  + k. 


4 H- 


(c2 


+ 


- 4c2t>jA2  ( k2x  + k2  + 


4 H 


2l  + 4g(g-  H)k 


2 

lJT 


A2  *F  Af  + 

l/2s 


(16) 


We  shall  refer  to  these  two  modes  as  the  ‘plus'  and  ‘minus'  modes.  These  modes 
involve  the  interaction  of  all  three  restoring  forces:  buoyant,  pressure,  and  magnetic. 
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They  reduce  to  more  familiar  wave  modes  in  certain  limiting  cases.  For  example,  in 
the  limit  of  vanishing  magnetic  field  (rA  — » 0)  they  reduce  to  the  acoustic  and  gravity 
modes  of  an  isothermal  atmosphere  (see,  for  example.  Thomas  a al.,  1971).  Alter- 
natively, in  the  limit  of  no  stratification  (//  -*  x , iy  -»0),  they  reduce  to  the  fast  and 
slow  magnetoacoustic  waves  in  a homogeneous  atmosphere  (see,  for  example,  Oster- 
brock,  1961).  For  intermediate  cases  (such  as  for  penumbral  conditions),  we  can  look 
at  the  modes  as  being  either  magnetically  modified  acoustic-gravity  waves  or  gravity- 
modified  magnetoacoustic  waves.  We  shall  continue  to  use  the  terms  plus  and  minus 
modes  here. 

The  plus  and  minus  modes,  for  k.=  0,  are  shown  schematically  in  a diagnostic 
diagram  in  Figure  t.  The  asymptotic  behavior  of  the  dispersion  relation  (16)  for 
k2  = 0 is  such  that 


Fig.  I.  Schematic  diagnostic  diagram  of  the  plus  and  minus  modes,  showing  curves  for  k;  -0 
(solid  lines).  The  plus  mode  has  a finite  cutoff  frequency  as  kx^> 0.  For  kx  -*■ c,  the  plus  and  minus 
modes  approach  the  dispersion  lines  for  the  pure  acoustic  and  pure  Alfven  modes.  For  strong  stratifi- 
cation and  weak  magnetic  field,  the  minus  mode  approaches  the  line  c>  - kxvx  from  above  rather 

than  belt  *. 
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In  the  diagnostic  diagram  (Figure  I),  there  are  vertically-propagating  waves  (A2  >0) 
above  the  curve  and  below  the  w_  curve,  and  evanescent  waves  (kl< 0)  in  the 
region  between  the  two  curves. 

The  plus  mode  has  a finite  cutoff  frequency  of  (oi+)(.=  (c2  + v\)'l2j2H  as  A,-»0. 
Estimates  of  this  cutoff  frequency  at  photospheric  levels  in  a penumbra  give  values 
just  a little  lower  than  the  frequencies  of  penumbral  waves.  This  suggests  that  the 
penumbral  waves  might  be  identified  with  plus  modes.  We  shall  show  this  in  more 
detail  through  the  use  of  our  penumbral  model. 


3.  The  Penumbral  Model 


We  now  present  the  penumbral  model  to  be  used  in  our  wave  calculations.  For  our 
purposes,  a penumbral  model  consists  of  specified  distributions  of  c2,  t*,  and  //  with 
height  z.  To  simplify  the  calculations,  we  have  chosen  to  represent  the  expected  verti- 
cal variation  of  these  parameters  by  piecewise  linear  functions.  The  specified  forms 
of  c2(z),  (- ).  and  //  (z)  are  show  n in  Figure  2.  The  important  features  of  the  model 

are  as  follows.  The  sound  speed  increases  with  depth  into  the  convection  zone  from 
a broad  minimum  in  the  penumbral  photosphere  and  chromosphere.  The  Alfven 
velocity  increases  rapidly  with  height  due  to  the  nearly  exponential  decrease  of 
density  while  the  magnetic  field  decreases  nearly  linearly.  The  density  scale  height 
increases  on  either  side  of  a minimum  in  the  penumbral  photosphere  and  low  chromo- 
sphere. The  model  thus  reproduces  the  main  expected  features  of  vertical  penumbral 
structure.  The  numerical  values  of  the  parameters  were  chosen  to  represent  a typical 
penumbra. 

The  distribution  of  the  sound  speed  was  determined  primarily  from  the  expected 
temperature  distribution.  Kjeldseth  Moe  and  Maltby  (1969)  report  that  the  tempera- 
ture T in  t.'-.c  penumbra  may  be  obtained  by  adding  a constant  A0  = 0.055  to  the  0 
values  of  the  quiet  photosphere,  where  0=  5040,  T.  With  this  A0,  the  calculated  relative 
intensities  averaged  over  the  penumbral  fine  structure  agreed  well  with  observations. 
Kjeldseth  Moe  and  Maltby  used  the  Bilderberg  Continuum  Atmosphere  (Gingcrich 
and  de  Jager,  1968)  for  their  quiet  photospheric  temperatures:  here,  we  use  the  more 
recent  Harvard-Smithsonian  Reference  Atmosphere  (HSRA,  Gingcrich  ct  al.,  1971). 
Using  the  constant  J0  = O.O55,  the  penumbral  temperature  minimum  is  found  to  be 
3989 K,  yielding  a minimum  value  of  the  sound  speed  squared  of  c2  = 43.5  '.m2  s'2. 
The  slope  of  the  linear  increase  of  c2  into  the  convection  zone  was  taken  to  be  0.1 
km  s"2;  this  choice  is  also  based  on  the  behavior  of  the  HSRA.  The  distribution  of 
sound  speed  (in  km  s"1)  is  given  by 


c2(z ) 


(63.5 -O.lz 
(43.5 


z 200 
z > 200. 


(19) 


At  some  point  in  the  penumbral  chromosphere  the  sound  speed  will  increase  rapidly 
to  coronal  values.  However,  this  increase  takes  place  above  the  region  where  the 
penumbral  waves  are  trapped  in  our  model,  so  we  have  not  included  this  in  the  model. 


11 


406 


ALAN  H.NYF  AND  JOHN  H. THOMAS 


z*  10'2  (km) 


z « 10' 2 (km) 

Fig.  2.  Distribution  of  c2,  rA2,  and  H with  height  : in  the  penumbral  model  (solid  lines).  The  data 
(crosses)  for  iA2  are  based  on  penumbral  observations  (see  text).  Data  points  for  quiet  photospheric 
values  of  r2  and  //,  based  on  the  HSRA,  are  shown  oniy  for  comparison  with  the  penumbral  model. 

The  distribution  of  c 2 is  shown  in  Figure  2,  along  with  quiet  atmospheric  values  of  c2 
from  the  HSRA  (Nakagawa,  1973)  for  comparison. 

To  determine  the  Alfven  velocity,  we  assumed  that  the  magnetic  field  strength 
decreases  linearly  with  height  and  is  1000  G at  z = 0.  The  rate  of  decrease  of  field 
strength  with  height  was  taken  to  be  0.2  G km-1,  in  accordance  with  observations  (see 
Bray  and  Loughhead,  1964).  The  photospheric  densities  were  estimated  from  the 
penumbral  model  of  Makita  (1963),  the  only  penumbral  model  to  give  densities.  His 
densities  are  consistent  with  more  recent  quiet  photospheric  models.  The  resulting 
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values  of  il  are  plotted  in  figure  2.  The  vertical  distribution  of  these  points  was 
approximated  in  our  model  by  two  straight-line  segments,  with  the  upper  segment 
extrapolated  to  greater  heights.  In  the  convection  zone,  the  Alfven  velocity  was 
approximated  by  another  straight  line  segment  which  goes  to  zero  at  z = -500  km, 
due  primarily  to  the  increase  in  density  with  depth,  but  also  to  the  fact  that  the  penunt- 
bral  field  probably  lies  over  the  convection  zone.  The  resulting  distribution  of  the 
Alfven  velocity  (in  km  s'1)  is  given  by 


i-a(-) 


0 

25  + 0.05z 

25  + O.I72z 

94  + 0.75  (z  - 400) 


z < - 500 
- 500  z < 0 
0 < z < 400 
z > 400. 


(20) 


The  distribution  of  the  local  density  scale  height  H (in  km)  with  height  was  chosen 
to  be 

( 230  - 0.375z  z < 200 

H (z)  = 155  200  ^ z < 400  (21) 

| 155  +0.55  (z- 400)  z>  400. 


In  choosing  values  for  H.  we  were  guided  by  the  following  expression  for  the  vertical 
entropy  gradient  in  the  atmosphere: 


d.'  = QC..K  r/  1>a\  1 

dz  PT  [V  2 / H 


(22) 


where  c„  is  the  specific  heat  at  constant  volume,  k is  the  isothermal  compressibility, 
and  p is  the  coefficient  of  thermal  expansion.  This  relation  follows  from  the  hydro- 
static equation  (1)  and  the  basic  thermodynamic  relation 


Tdi  = 


KC, 

P 


(dp  — c"  dp). 


(23) 


The  distribution  of  H given  in  (21 ) is  such  that  the  entropy  gradient  dv  dz  is  negative  in 
the  upper  convection  zone  and  positive  in  the  penumbral  photosphere  and  chromo- 
sphere. The  distribution  of  H is  shown  in  Figure  2,  along  with  quiet  atmospheric 
values  of  H from  the  HSRA  (Nakagawa,  1973)  for  comparison.  Here  we  can  see  the 
effect  of  the  supporting  magnetic  field;  the  scale  height  is  greater  in  the  penumnra 
than  in  the  quiet  atmosphere.  The  effect  of  the  magnetic  field  increases  as  the  density 
decreases. 

4.  Trapping  of  Plus  Modes  in  the  Penumbral  Phoiosphere 


We  now  show  that  the  running  penumbral  waves  can  be  identified  in  our  model  with 
plus  modes  that  are  vertically  trapped  at  photosphenc  levels.  The  basic  mechanism  for 
the  trapping  is  the  refraction  due  to  (i)  the  increasing  Alfven  velocity  with  increasing 
height  in  the  photosphere-low  chromosphere,  and  (ii)  the  increasing  sound  speed  with 
depth  in  the  convection  zone.  V/e  shall  demonstrate  the  trapping  mechanism  in  two 
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ways:  first,  by  considering  the  local  dispersion  relation  (16)  at  three  distinct  levels  in 
the  model  penumbra  (assuming  locally  constant  parameters),  and  second,  by  con- 
sidering the  propagation  equation  (10)  for  the  complete  continuous  penumbral 
model  (Figure  2). 

4.1.  Three-level  model 

A convenient  method  of  looking  at  wave  modes  in  our  model  is  to  draw  diagnostic 
diagrams  such  as  Figure  1 for  various  heights  in  the  atmosphere,  assuming  locally 
constant  values  of  c 2,  ti,  and  H - that  is,  ignoring  the  derivatives  of  these  parameters 
in  Equation  (3)— (5).  Although  only  approximate,  this  method  does  give  some  feeling 
for  the  behavior  of  the  wave  modes.  We  shall  draw  diagnostic  diagrams  for  three 
different  levels  in  our  model  penumbra;  one  in  the  convection  zone,  one  in  the 
photosphere,  and  one  in  the  chromosphere.  We  may  interpret  these  diagnostic 


\ « 10  3 (krn) 

5 3 2 1.5 


kx » 10 J (km-1 ) 

Fig.  3.  Superimposed  diagnostic  diagrams  for  values  of  r 2,  r*1.  and  H at  three  different  levels  in 
the  model  penumbra:  level,  I z~  - 1400  km;  level  2,  : 50  km;  level  3.  : 500  km.  The  shaded 

region  is  a region  of  vertical  trapping  of  plus  modes  around  level  2.  The  crosses  correspond  to 
observed  running  penumbral  waves  (Giovanelli,  1974). 
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diagrams  more  properly  as  representing  a three-layer  model  of  the  penumbra,  where 
in  each  layer  the  parameters  are  constant.  A middle  layer  of  finite  thickness  represent- 
ing the  penumbral  photosphere  is  bounded  above  and  below  by  semi-infinite  layers 
representing  the  penumbral  chromosphere  and  convection  /one.  respectively. 

Figure  3 shows  diagnostic  diagrams  for  values  of  c2,  rj,  and  //  at  three  different 
levels  in  the  penumbral  model:  1400  km  (layer  1.  convection  zone).  z = -50  km 

(layer  2.  penumbral  photosphere),  and  z = 50()  km  (layer  3,  penumbral  chromosphere). 
The  important  feature  of  this  figure  is  the  existence  of  the  shaded  region  in  which  the 
plus  modes  are  vertically  propagating  in  layer  2,  but  are  vertically  evanescent  in 
layers  I and  3.  These  modes  are  thus  trapped  in  layer  2.  The  plus  mode  diagnostic- 
curve  for  layer  3 lies  above  that  of  layer  2 because  of  the  higher  Alfven  velocity  in 
layer  3.  For  higher  levels  in  the  chromosphere,  the  plus  mode  curve  for  layer  3 will  be 
higher  than  that  shown.  The  plus  mode  diagnostic  curve  for  layer  I has  a lower  cutoff 
frequency  than  layer  2 due  to  the  larger  scale  height  in  the  convection  zone,  but  lies 
above  the  curve  for  layer  2 for  higher  values  of  A , because  of  the  higher  value  of  c2  in 
layer  I . For  deeper  levels  in  the  convection  zone,  the  plus  mode  curve  for  layer  I will 
have  an  even  lower  cutoff  frequency  and  a steeper  slope.  The  existence  of  the  region  of 
trapping  in  the  diagnostic  diagram  is  a consequence  of  the  qualitative  features  of  our 
model  and  is  not  dependent  on  the  particular  choice  of  numerical  values. 

There  is  no  trapping  of  minus  modes  in  Figure  3.  All  the  minus  modes  propagating 
in  layer  2 are  also  propagating  in  layer  3.  The  minus  mode  curve  for  layer  1 has 
vanished  since  the  Alfven  velocity  is  zero  and  the  atmosphere  is  convectivelv  unstable 
at  that  level. 

Also  shown  in  Figure  3 are  four  data  points  corresponding  to  penumbral  waves 
observed  in  different  sunspots  by  Giovanelli  (1974)  for  which  he  gives  specific  wave- 
length and  periods.  These  points  tend  to  cluster  in  the  long-wavelength  end  of  the 
shaded  region  of  trapping  of  the  plus  modes.  The  one  point  which  lies  outside  the 
shaded  region  does  lie  in  the  region  of  trapping  if  the  parameter  values  in  layer  I are 
chosen  to  represent  a lower  level  in  the  convection  zone.  That  is.  a wave  of  this  fre- 
quency and  wavelength  is  reflected  at  a lower  level  in  the  convection  zone  in  our 
model.  This  data  point  is  considered  atypical  by  Giovanelli,  however.  The  observation- 
al data  correspond  to  dilferent  sunspots,  no  doubt  having  different  field  strengths, 
whereas  the  diagnostic  curves  are  for  a single  choice  of  the  model  parameters.  Never- 
theless, the  resulting  picture  in  Figure  3 clearly  shows  that  the  running  penumbral 
waves  should  be  identified  with  plus  modes  which  are  vertically  trapped  at  photo- 
spheric  levels  in  our  model. 

4.2.  Continuous  model 

We  now  illustrate  more  accurately  the  trapping  of  the  penumbral  waves  by  making  use 
of  our  complete  penumbral  model  (Figure  2)  in  conjunction  with  the  propagation 
equation  (10).  Roughly  speaking,  a wave  is  vertically  propagating  when  the  function 
/ (At,  to:  z)  (Equation  ( 1 1 ))  is  positive,  and  is  vertically  evanescent  when /is  negative. 
In  Figure  4 we  have  plotted /as  a function  of  z for  our  penumbral  model  (Figure  2)  for 
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z ■ 10’*  ( km ) 


Fig.  4.  The  function /(A,,  oi;  :)  evaluated  for  the  penumbral  model  (Figure  2)  with  kx  2 10  3 

km"  cj  3 x 10”  1 s~ 


a typical  value  of  horizontal  wavenumber  (kx  = 2 x 10' 3 km ) and  frequency 
(co  = 3x  10~ 2 s_  1 ) for  penumbral  waves.  The  function  f (z)  is  discontinuous  due  to 
the  discontinuities  in  the  first  derivatives  of  the  parameters  c2,  r*.  and  H in  our  model. 
The  general  behavior  of  / is  that  it  is  positive  in  a central  region  extending  from 
z*  — 1000  km  to  roughly  z%300  km,  and  negative  above  and  below  this  region. 

In  the  lowest  order  WKB  approximation,  the  general  solution  to  ( 10)  is  given  by 


4>(z)  = exp 


(24) 


If /is  positive  for  0<z<z0  and  negative  for  z>z0,  then  the  solution  for  z>z0  which  is 
bounded  as  z->  oo  can  be  written  as 


<t>  (z)  = exp£±  i | Jfdz-  J^-ZdzJ.  (25) 

0 ro 

Thus,  for  z>z0,  the  solution  has  the  form  of  an  exponential  decay  (with  variable 
exponential  factor),  and  the  wave  may  be  characterized  as  evanescent  in  the  vertical 
direction.  This  analysis  is  somewhat  crude,  however,  since  the  WKB  solution  (24)  is 
not  accurate  in  the  neighborhood  of  z0.  The  same  analysis  applies  to  the  case  where/ 
is  positive  for  — Zq<z<0  and  negative  for  z<  —z„.  The  wave  then  becomes  evanes- 
cent as  z decreases  below  — z'0.  Applying  this  to  Figure  4.  we  conclude  that  a wave 
with  the  specified  horizontal  wavelength  and  frequency  is  trapped  roughly  in  the 
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central  region  — 1000<z<300  where /is  positive.  The  exact  location  of  the  upper 
boundary  of  the  trapping  region  is  somewhat  unclear  due  to  the  complicated  dis- 
continuous behavior  of/ there.  This  is  due  to  the  piecewise  linear  nature  of  our  model; 
with  smoother  distributions  of  c2,  v and  //,  the  function  / would  have  a single, 
smooth  zero  crossing  at  some  z>0. 

We  have  been  discussing  the  behavior  of  the  function  /(z).  Comparison  with  the 
case  of  constant  parameters  shows  that  |0|2  is  roughly  proportional  to  the  kinetic- 
energy  of  vertical  motion  of  the  wave  (<p  ~ sJq  iv).  Thus,  it  is  the  kinetic  energy  that  is 
trapped  in  the  region  of  positive  / The  amplitude  of  the  vertical  velocity  behaves 
somewhat  differently.  From  Equations  (9)  and  (25)  we  have,  for  z>z0. 


2° 

w 


li-(z)  = exp  ± i ^//dz  + Q(z) 


z z 

Q(z)  -i  J/ldz-  JV-Zdz. 


The  function  A (z)  (Equation  (7))  is  negative  over  a range  in  z,  and  the  function  Q (z) 
will  remain  positive  over  a range  of  z beyond  z0  before  becoming  negative.  Thus,  the 
vertical  velocity  continues  to  increase  above  z0  before  reaching  a maximum  at  z = z() 
where  z,  is  determined  from  the  relation  Q(zl)  = 0.  This  can  occur  in  spite  of  the 
decreasing  energy  above  z = z0  due  to  the  rapid  decrease  in  density. 

By  numerically  evaluating  the  integrals  in  (27),  we  have  estimated  the  point  z,  of 
maximum  vertical  velocity  in  our  model,  again  for  the  choice  kx  = 2 x 10” J km-1  and 
w = 3 x 1(T2  s~‘.  The  resulting  value  is  z,  s:  1250  km.  This  is  at  a level  where  the  wave 
is  evanescent,  and  is  well  into  the  region  of  formation  of  Fla  (Vernazza  et  al.,  1973). 
Thus,  although  the  trapped  waves  have  their  maximum  energy  at  lower  levels  (z%0) 
where  the  density  is  higher,  the  vertical  velocity  is  greatest  at  chromospheric  levels 
where  is  it  observed. 

5.  Discussion 

An  important  conclusion  that  emerges  from  our  model  is  that  the  running  penumbral 
waves  are  basically  a photospheric  phenomenon,  even  though  they  are  observed  at 
chromospheric  levels.  This  is  really  an  expected  result;  the  phase  velocity  of  the 
penumbral  waves  is  typical  of  photospheric  conditions  rather  than  chromospheric. 
For  example,  the  Alfven  velocity  at  Fiat  levels  in  the  penumbra  is  far  greater  than  the 
phase  velocity  of  the  penumbral  waves.  In  our  picture,  the  waves  are  evanescent  where 
they  are  observed.  That  is,  the  observed  vertical  motions  are  a passive  response  to  an 
actively  propagating  wave  at  lower  levels.  We  have  shown  that  although  the  energy  is 
trapped  at  photospheric  levels,  the  vertical  velocity  reaches  its  maximum  value  at 
levels  of  formation  of  Fix,  consistent  with  observations. 

The  concept  of  the  penumbral  waves  as  a photospheric  phenomenon  is  also  quite 
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consistent  with  their  expected  source  of  excitation.  Moore  (1973)  has  shown  that  the 
likely  source  of  excitation  of  the  penumbral  waves  is  oscillatory  convection  in  a sub- 
photospheric  layer  in  the  umbra.  The  observation  that  the  penumbral  waves  dis- 
appear suddenly  at  the  boundary  between  the  white-light  penumbra  and  the  sur- 
rounding photosphere  also  supports  the  photospherie  nature  of  the  waves.  In  the 
chromosphere,  the  penumbral  libril  structure  extends  outward  beyond  this  boundary. 

The  basic  mechanism  for  the  vertical  trapping  of  the  penumbral  waves  is  not  sensi- 
tive to  the  details  of  our  penumbral  model.  This  mechanism  is  based  primarily  on  the 
increase  in  the  sound  speed  as  we  go  down  into  the  convection  /one  and  the  increase 
in  the  Alfven  velocity  as  we  go  up  into  the  chromosphere. These  features  are  certain 
to  remain  in  any  improved  penumbral  model.  Observationally.  the  insensitivity  of  the 
basic  mechanism  is  confirmed  by  the  fact  that  penumbral  waves  are  seen  in  almost 
every  stable,  regular  sunspot  (Zirin  and  Stein.  1972). 

The  character  of  the  trapped  plus  modes  is  different  at  different  heights.  The  wave  is 
more  nearly  acoustic  at  low  levels  (convection  /one),  but  is  more  nearly  Alfvenic  at 
higher  levels  (photosphere  and  low  chromosphere).  The  effect  of  stratification  and 
gravity,  while  not  dominant,  is  felt  throughout  the  trapping  region. 

Although  we  have  not  attempted  to  calculate  them,  there  are  certainly  resonant 
modes  in  our  model  which  arise  from  the  constructive  interference  of  plus  modes 
reflected  from  above  and  below  in  the  trapping  region.  The  resonant  dispersion  curves 
would  lie  in  the  shaded  region  in  Figure  3.  However,  since  the  penumbral  waves  are  a 
somewhat  transitory  phenomenon,  there  is  no  reason  to  associate  them  specifically 
with  such  a resonant  mode:  it  is  sufficient  that  a wave  of  the  proper  frequency  and 
horizontal  wavelength  be  vertically  trapped. 

Since  our  penumbral  model  is  horizontally  uniform,  it  cannot  account  for  any 
observed  horizontal  (radial)  variations  in  the  penumbral  waves.  We  can  roughly 
account  for  the  radial  geometry  of  the  penumbra  by  noting  that  the  energy  in  a verti- 
cally trapped  wave  will  decrease  as  I V.  where  r is  radius  measured  from  the  spot 
center,  as  the  wave  propagates  radially  outward.  This  accounts  at  least  in  part  for  the 
observed  decrease  in  wave  amplitude  with  radius.  The  sudden  disappearance  of  the 
waves  at  the  outer  edge  of  the  penumbra  is  associated  with  the  more  fundamental 
problem  of  the  basic  existence  of  the  sharp  boundary  between  the  white-light  penum- 
bra and  the  surrounding  photosphere. 

Finally,  we  should  note  that  our  assumption  of  adiabatic  perturbations  is  invalid 
over  a limited  range  of  height  in  the  low  photosphere,  say  the  first  few  hundred 
kilometers  above  c = 0.  Here  one  should  account  for  the  rapid  radiative  exchange.  The 
main  effect  of  the  radiative  transfer  on  the  plus  modes  is  to  effectively  replace  the 
adiabatic  sound  speed  in  this  limited  region  by  something  nearer  to  the  isothermal 
sound  speed.  The  overall  effect  on  the  trapped  plus  modes  discussed  here  will  be  quite 
small,  even  quantitatively. 
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Convective  instability  in  the  presence  of  a nonuniform 
horizontal  magnetic  field 
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Newcomb’s  criterion  for  convective  stability  in  the  presence  of  a horizontal  magnetic  field  is  written 
in  a form  which  explicitly  shows  the  effect  of  vertical  variations  of  the  magnetic  field  strength.  It  is 
shown  that  a nonuniform  horizontal  magnetic  field  can  be  destabilizing  as  well  as  stabilizing. 


We  consider  the  convective  instability  of  a compres- 
sible, inviscid,  perfectly  conducting  gas  permeated  by 
a horizontal  magnetic  field  B = [Mz),  0, 0]  that  may  vary 
with  height  z,  under  a uniform  gravitational  accelera- 
tion g (in  the  negative  z direction),  assuming  adiabatic 
perturbations.  The  first  complete  treatment  of  this 
problem  was  given  by  Newcomb1  using  the  energy  inte- 
gral method.  He  showed  that  a necessary  and  sufficient 
condition  for  stability  is  given  by 


>P^ 

dz  yp  ’ 


(1) 


where  p is  mass  density,  p is  pressure,  and  y is  the 
ratio  of  specific  heats.  In  the  case  of  instability,  the 
most  unstable  mode  has  the  form  of  an  interchange  of 
long  but  finite  segments  of  magnetic  field  lines.  The 
stability  criterion  (1)  was  also  derived  by  Yu2  by  con- 
sidering the  force  balance  on  a displaced  magnetic  flux 
tube. 

Newcomb  noted  that  the  critical  density  gradient  on 
the  right-hand  side  of  (1)  is,  at  least  explicitly,  inde- 
pendent of  the  magnetic  field,  and  is,  in  fact,  that  given 
by  the  Schwarzschild  criterion3  in  the  absence  of  a mag- 
netic field.  However,  in  the  case  of  a nonuniform  hori- 
zontal magnetic  field  (for  which  Newcomb’s  analysis  is 
valid),  the  static  distribution  of  pressure  and  density  is 
affected  by  the  magnetic  field,  and  thus  the  stability 
criterion  (1)  depends  implicitly  on  the  magnetic  field. 
This  point,  which  was  overlooked  by  Newcomb,  will  be 
pursued  here. 


Using  the  equation  of  hydrostatic  equilibrium, 


and  the  equation  of  state  of  a perfect  gas,  p = pRT,  the 
stability  criterion  (1)  can  be  written  in  the  form 


dT_  /<m  > _ _1_ 

dz  \rfz  /,  pR  dz  \8t  / ’ 


(3) 


where  (d T/dz),  = - g/c,  = - g(y  - \)/>R  is  the  adiabatic 
temperature  gradient.  The  new  form  (3)  of  the  stability 
criterion  has  an  advantage  over  the  form  (1)  in  that  it 
shows  explicitly  the  effect  of  vertical  variations  of  the 
magnetic  field,  since  the  temperature  distribution  is 
independent  of  the  magnetic  field.  For  a uniform  mag- 
netic field  ( B=  const),  (3)  reduces  to 


which  is  identical  to  the  Schwarzschild  criterion  in  the 
absence  of  a magnetic  field.  [Note  that  usually,  e.g. , 
in  a star,  dT/dz  and  ( dT/dz)s  are  both  negative.]  Thus, 
a uniform  horizontal  magnetic  field  has  no  effect  on  the 
condition  for  the  onset  of  convective  instability, 4 al- 
though, as  Newcomb1  also  showed,  it  does  have  an  ef- 
fect on  the  growth  rates  of  unstable  modes. 

If,  however,  the  magnetic  field  is  nonuniform,  then 
(3)  shows  that  the  field  can  be  stabilizing  (in  the  case 
dB/dz  >0)  or  destabilizing  ( dii/dz  <0).  A field  that  in- 
creases with  height  (dB/dz  > 0)  can  stabilize  the  atmo- 
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sphere  in  the  presence  of  a superadiabatic  temperature 
gradient.  When  the  field  decreases  with  height  ( dB/dz 
<0),  the  critical  temperature  gradient  is  reduced  be- 
low the  adiabatic  gradient. 

As  an  illustration  of  the  destabilizing  effect  for  dB / 
dz  <0,  consider  the  case  of  an  isothermal  atmosphere, 
which  is  convectively  stable  in  the  absence  of  a mag- 
netic field.  The  stability  criterion  (3)  reduces  to 


<_*(*!) -Jyjil) . 

p dz  \8tt / \dz)s  K\y) 


(4) 


If  the  magnetic  pressure  decreases  more  rapidly  with 
height  than  the  critical  rate 


[=©L— Or1)- 


(5) 


the  atmosphere  will  be  convectively  unstable.  For  this 
critical  magnetic  field  gradient,  the  corresponding 
pressure  gradient  is  found  from  (2)  to  be 


\dj>]  _ _£g_ 

LdzJc  y 


(6) 


Thus,  in  the  state  of  marginal  stability  in  an  isothermal 
atmosphere,  the  gradient  of  gas  pressure  balances  the 
fraction  l/y  of  the  gravitational  force,  while  the  grad- 
ient of  magnetic  pressure  balances  the  remaining  frac- 
tion (y  - l)/y  of  the  gravitational  force. 

A special  case  of  interest  is  that  of  an  isothermal  at- 
mosphere in  which  the  field  decreases  with  height  in 
just  the  manner  that  causes  the  Alfven  speed  va  = {B?/ 
4irp)l/2  to  remain  constant.  In  this  case,  (4)  can  be 
written  as 


^ ) ■ 


(7) 


where  c=(yRT)l,z  is  the  adiabatic  sound  speed.  This 
stability  criterion  was  obtained  for  this  special  case  by 
Yu5  from  a normal  mode  analysis.  Parker6  has  studied 
this  case,  including  the  effect  of  cosmic -ray  pressure, 
in  connection  with  the  gaseous  disk  of  the  galaxy. 

The  destabilizing  effect  of  a nonuniform  magnetic 
field  with  dB/dz  < 0 is  similar  to  the  phenomenon  of 
“magnetic  buoyancy”  analyzed  by  Parker7  and  by  Weiss.8 
Magnetic  buoyancy  is  attributed  to  an  isolated  tube  of 
magnetic  flux  in  thermal  equilibrium  with  its  nonmag- 
netic surroundings.  The  instability  discussed  here  is 
due  to  the  buoyancy  of  an  arbitrary  tube  of  flux  within 
a smoothly  varying  magnetic  field. 
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Predoctoral  Fellowship. 
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ABSTRACT 


The  linearized  theory  of  magneto-atmospheric  waves  (involving 
the  combined  restoring  forces  due  to  buoyancy,  compressibility, 
and  magnetic  field)  is  developed  for  the  case  of  a horizontal 
magnetic  field.  A general  propagation  equation  is  derived  for 
adiabatic  perturbations  with  arbitrary  vertical  distributions  of 
the  sound  speed  c,  Alfven  velocity  v^,  and  local  density  scale 
height  H.  An  exact  analytical  solution  to  the  propagation 
equation  is  obtained  for  the  case  of  an  isothermal  atmosphere 
permeated  by  a uniform  horizontal  magnetic  field,  without  making 
the  usual  short-wavelength  assumption.  This  solution  is  applied 
to  an  idealized  model  of  the  low-corona-chromosphere  transition 
region  for  comparison  with  observations  of  flare- induced  coronal 
waves.  The  results  show  that  disturbances  may  propagate  horizon- 
tally in  the  low  corona  in  a wave  guide  formed  by  the  sudden 
density  increase  into  the  chromosphere  below  and  by  the  rapidly 
increasing  Alfven  velocity  with  height  in  the  corona.  The  group 
velocity  of  the  guided  wave  modes  is  nearly  independent  of  wave- 
length, so  that  a disturbance  propagates  as  a compact  wave  packet. 
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I.  INTRODUCTION 

The  theory  of  waves  in  a compressible,  stratified, 
electrically-conducting  atmosphere  permeated  by  a magnetic 
field  is  of  considerable  importance  in  astrophysics,  especially 
in  solar  physics  where  there  is  a wealth  of  detailed  observations 
of  such  waves  in  the  solar  atmosphere.  Following  Yu  (1965),  we 
shall  refer  to  waves  which  involve  compressibility,  buoyancy, 
and  magnetic  forces  as  magneto-atmospheric  waves.  Among  the 
many  solar  phenomena  that  are  seemingly  attributable  to  magneto- 
atmospheric  waves  are  the  heating  of  the  corona  and  of  chromo- 
spheric plages,  the  five-minute  oscillations  in  active  regions, 
oscillations  in  sunspot  umbrae  and  penumbrae,  and  flare-induced 
coronal  disturbances. 

The  theory  of  magneto-atmospheric  waves  is  complicated  by  the 
anisotropic  nature  of  the  medium;  the  gravitational  field  and 
the  magnetic  field  each  introduce  preferred  directions.  Addition- 
ally, the  disturbance  is  subjected  to  the  combined  restoring 
forces  due  to  compressibility,  buoyancy,  and  the  magnetic  field, 
so  that  pure  wave  modes  (i.e.,  acoustic,  Alfven,  and  gravity) 
exist  only  as  special  cases.  In  general,  magneto-atmospheric 
waves  involve  the  effects  of  all  three  restoring  forces.  The 
problem  is  compounded  by  the  fact  that  the  basic  parameters 
describing  wave  propagation  in  the  solar  atmosphere  (the  sound 
speed  c,  the  Alfven  velocity  v^,  and  the  local  density  scale  height 
H)  are,  in  general,  functions  of  height,  and  therefore  the 
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disturbance  cannot  be  represented  by  plane  waves  propagating 
in  the  atmosphere. 

Previously,  the  problem  of  magneto-atmospheric  waves  has 
been  studied  by  either  of  two  basic  approaches,  each  yielding 
a dispersion  relation  based  on  constant  values  of  the  atmospheric 
parameters.  The  first  approach  has  been  to  assume  that  the 
vertical  extent  of  the  disturbance  is  much  less  than  the  smallest 
scale  height  for  variation  of  the  atmospheric  parameters.  The 
parameters  can  then  be  taken  as  constant  locally.  McLellan  and 
Winterberg  (1968)  studied  an  isothermal  atmosphere  permeated  by 
a uniform  magnetic  field  with  arbitrary  orientation.  Then, 
assuming  that  the  Alfven  velocity  is  constant  locally  (although 
it  actually  increases  exponentially  with  height) , they  derived 
a local  dispersion  relation  that  is  valid  for  short  wavelengths 
(short  compared  to  the  density  scale  height)  . This  local  dis- 
persion relation  has  been  studied  by  several  authors  (Bel  and 
Mein  1971;  Michalitsanos  1973;  Nakagawa,  Priest,  and  Wellck 
1973)  to  determine  the  effects  of  different  propagation  direc- 
tions and  magnetic  field  orientations. 

The  second  basic  approach  has  been  to  investigate  an  iso- 
thermal atmosphere  permeated  by  a horizontal  magnetic  field  that 
decreases  exponentially  with  height  such  that  the  Alfven  velocity 
is  constant,  as  are  the  sound  speed  and  the  density  scale  height. 
Yu  (1965)  derived  the  dispersion  relation  which  in  this  case  is 
valid  for  all  wavelengths.  He  evaluated  the  three  modes  of 
propagation  for  various  angles  between  the  wave  propagation  vector 


and  the  magnetic  field.  Chen  and  Lykoudis  (1972)  used  the 
dispersion  relation  of  Yu  to  study  the  five-minute  oscilla- 
tions in  plage  regions.  Nye  and  Thomas  (1974a)  used  Yu's 
dispersion  relation  in  connection  with  a multi-layer  model 
for  running  penumbral  waves. 

We  consider  only  plane  parallel  atmospheres,  with  no 
horizontal  variation,  in  order  to  permit  Fourier  decomposition 
in  the  horizontal  directions  as  well  as  time.  With  variations 
in  the  z-direction  only,  the  governing  partial  differential 
equations  reduce  to  ordinary  differential  equations.  There 
are  only  two  magnetic  field  configurations  consistent  with 
static  equilibrium  and  no  horizontal  variation.  They  are 
B = constant  and  B = (B^ (z) ,B  (z)  , 0)  . We  shall  restrict  our 
study  to  the  case  of  a unidirectional  horizontal  magnetic  field 
that  may  vary  with  height.  We  shall  not  follow  either  of  the 
two  basic  approaches  discussed  above,  however,  due  to  the 
inherent  limitations  of  each  of  them  for  solar  applications. 

The  length  scale  for  magnetic  field  changes  in  the  solar  atmo- 
sphere is  generally  much  greater  than  the  density  scale  height, 
so  that  the  assumption  of  constant  Alfv^n  velocity  is  not 
justified.  On  the  other  hand,  the  length  scale  of  observed 
disturbances  in  the  solar  atmosphere  is  not  generally  small  com- 
pared to  the  density  scale  height.  This  is  especially  true  in 
the  photosphere  where  the  density  scale  height  may  be  smaller 
than  100  km,  and  hence  smaller  than  the  limit  of  observational 


resolution. 


In  section  II  we  derive  a general  propagation  equation 
for  an  arbitrary  direction  of  propagation  and  arbitrary 
vertical  variations  of  the  atmospheric  parameters.  For  the 
case  of  an  isothermal  atmosphere  permeated  by  a uniform  hori- 
zontal magnetic  field,  the  Alfvdn  velocity  increases  exponen- 
tially with  height.  In  section  III  we  obtain  an  exact  general 
solution  of  the  propagation  equation  in  this  case.  We  compute 
eigenmodes  for  the  case  of  a rigid  lower  boundary  in  section  IV, 
and  apply  this  to  a specific  solar  wave  phenomenon,  the  flare- 
induced  coronal  waves,  in  section  V. 

The  analytical  treatment  in  sections  II  and  III  also  forms 
the  basis  of  a following  paper  in  which  we  deal  with  running 
penumbral  waves. 


II.  BASIC  EQUATIONS 

The  atmosphere  is  assumed  to  be  a compressible,  inviscid, 
perfectly  conducting  gas  under  a uniform  acceleration  of 
gravity  g (=0.274  km  sec-2)  in  the  negative  z-direction.  The 
undisturbed  magnetic  field  is  taken  in  the  x-direction  and  may 
vary  with  height z;  i.e.,  B0  = (B0(z),O,O).  The  undisturbed 
pressure,  density,  and  temperature  may  all  be  functions  of 
height  z,  and  are  denoted  by  p0(z),  p0(z),  and  T0(z),  respectively. 
We  shall  see  that  wave  propagation  in  the  basic  atmosphere  may  be 
completely  characterized  by  the  vertical  variation  of  the  sound 
speed  c(z),  the  Alfven  velocity  vA(z),  and  the  local  density 


scale  height  H(z),  defined  by 


c 


2 


,3P< 
' 3p  ( 


Bo2  1 
4rrp0  ' H 


_1_  9P0 

p o 3 z 


(1) 


The  unperturbed  atmosphere  is  taken  to  be  in  static 
equilibrium: 
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If  the  magnetic  field  is  a function  of  height  z,  then  it  has  a 
role  in  the  basic  equilibrium  of  the  atmosphere.  We  consider 
only  stable  atmospheres,  which  requires  that 
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where  (dT0/dz)s  is  the  adiabatic  temperature  gradient  and  R is 
the  gas  constant  (see  Thomas  and  Nye  (1975)  for  a recent 
discussion) . 

Consider  small  adiabatic  perturbations  of  the  equilibrium 
atmosphere,  letting  p,  p,  u,  and  B denote  the  perturbations  in 
pressure,  density,  velocity,  and  magnetic  field,  respectively. 
Then,  the  basic  linearized  equations  of  continuity,  momentum, 
energy,  and  induction  are 

+ V • ( p ou)  = 0 , (3) 


Pol?  + Vp  “ p2.“  J ( VxB  o ) xB  + (VxB)xB0 


= 0, 
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3B 

- Vx (uxBfl ) = 0 . 


(6) 


After  taking  the  time  derivative  of  the  momentum  equation 
(4),  we  may  eliminate  the  perturbation  quantities  p,  p,  and  B 
by  using  equations  (3),  (5),  and  (6).  This  leaves  a single 

vector  equation  for  the  velocity  perturbation  u = (u,v,w): 


3 2u  1 i dB „ 2 , , I 

It?  + v(poW[g  + 8^7  "diT1  - c + 

+ J (u • V ) p o + PoV'uJg  - — j(VxBo) x [Vx (uxB0) ] - 
- B 0x (Vx  [Vx (uxB o ) ] ) | = 0. 


(7) 


Next,  we  assume  that  the  perturbation  velocity  has  the  form 
u = u exp  i(k«r-wt),  with  u = u(z)  = (u  (z)  , v (z)  , w (z)  ) and 
k*r  = k x + k y.  Then,  using  the  definition  of  va 2 , the  three 
components  of  the  momentum  equation  (7)  become 

(u)2-c2kx2)u  - c2kxkyv  - ikx  (g-c2^)  w = 0 , (8) 


- c2kxkyu  + [w2-c2ky2-vA2 (kx2+ky2) ] v - iky [g  - 

- <c!+V>a7!”  - 0 ■ 


and 


ik 


| dc ' 


dz 


- — + a + o2—  'u 

H q C dz  |u 


ik  |^-(c2 
y Idz 


+ v 2)  - 


(C2  + V.  ) 


A ' , _ . ,_2  . „ 2.  d I-  . ,_2  . - 2,  d2W 


+ q + (c  + VA  + (c  +VA 


* + - (c  + o^-v^ios  = o 


H I dz 


A x 


(9) 


(10) 


29 


They  give 


! 


These  are  the  linearized  perturbation  equations, 
important  information  about  particle  motions  for  various 
modes  of  propagation  in  the  atmosphere. 

The  horizontal  components  of  the  perturbation  velocity 
can  be  eliminated  from  the  system  of  equations  (8)  - (10)  to 


yield  a single  equation  for  the  vertical  velocity  w.  The 
resulting  equation  is 


d 2w 
dz  2 


A(z)  a? 


+ B (z)  w = 0 


(11) 


where  the  coefficients  A(z)  and  B(z)  are  given  by 
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Here,  D and  E are  given  by 


D = (to2-vA2kx2)[to2  (c2+vA2)-c2vA2kx2] 


(14) 


and 


E = u)4-(kx2+ky2)  [to2  (c2+vA2)-c2vA2kx2 


(15) 
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Equation  (11)  is  a general  equation  for  the  vertical  component 

of  velocity  for  a perturbation  propagating  in  an  arbitrary 

direction,  in  an  atmosphere  with  a horizontal  magnetic  field 

and  arbitrary  vertical  distributions  of  c2,  v^2,  and  II,  This 

propagation  equation  has  been  given  previously  (Nye  and  Thomas 

1974a)  in  the  case  k = 0. 

y 

III.  ISOTHERMAL  ATMOSPHERE  WITH  A UNIFORM 
HORIZONTAL  MAGNETIC  FIELD 


Now  consider  the  case  where  the  undisturbed  temperature 
and  magnetic  field  are  constant  with  height.  Since  the  magnetic 
field  is  uniform,  it  has  no  effect  on  the  hydrostatic  equilibrium 
of  the  atmosphere,  and  the  equilibrium  pressure  and  density  both 
decrease  exponentially  with  height.  The  sound  speed  and  the 
density  scale  height  are  both  constant,  with  values  determined 
by  the  temperature  of  the  atmosphere.  The  Alfv£n  velocity 
increases  exponentially  with  height  due  to  the  decreasing 
density.  The  sound  speed,  density  scale  height,  density,  and 
Alfv£n  velocity  are  given  by 

j. 

c = (yRTo) 2 = const.  , (16) 

RT. 

H = = const.  , (17) 

p0(z)  = Pooe“z/H  , (18) 

and 
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z/2H 


(19) 


vA ( z ) = v o e 


where  poo  and  v0  are  the  values  of  the  undisturbed  density  and 
Alfven  velocity  at  z = 0. 

The  nondimensional  parameter  6 2 = vo2/c2  is  introduced  as 
a measure  of  the  relative  importance  of  the  restoring  forces 
due  to  the  magnetic  field  and  to  compressibility  at  the  point 
z = 0.  For  values  of  (32  < 1,  there  is  a region  of  the  atmosphere 
above  z = 0 where  compressibility  has  more  importance  as  a 
restoring  force  than  does  the  magnetic  force,  but  in  any  case 
the  magnetic  field  always  becomes  dominant  as  z becomes  large. 

In  the  remainder  of  this  paper  we  shall  consider  only  waves 
whose  horizontal  component  of  propagation  is  parallel  to  the 
magnetic  field  (k  = 0).  Using  (16)  - (19),  we  may  write  (11)  as 


fc 


2w2  + 


u>2-c2k^)v0 


2 .z/H i d 2\ 

e Jdi" 


2 2 
C O) 


dw 

dz 


+ [(oi  -c 


2.o2i,2>  /,.2_,,2  Z/H.  2 , „ _ C,  . 2 


k^)  (ui2-v; 


kx)-q(g  -ir)kxlw  = 


(20) 


We  define  the  nondimensional  frequency  f)  and  the  non-dimensional 
horizontal  wavenumber  K by 


»■¥  , 


K = Hk 


(2i: 


By  transforming  the  dependent  and  independent  variables  according 
to 


W = we 


ZK/H 


x = 


62 (K2-fi2)  e 


-z/H 
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we  may  put  equation  (20)  in  the  dimensionless  form 
X(l-X)^£  + [C- (A+B+l)x]^  - ABW  = 0 , 


with 


A + B = C = 2K+  1 


dx 


AB  = ft2  + K + (*pr)|s 


(23) 


(24) 


Equation  (23)  is  the  standard  form  of  the  hypergeometric 
differential  equation.'*'  The  solutions  of  this  equation  may  be 
expressed  in  terms  of  hypergeometric  functions,  given  for 
| x | < 1 by 

n 


(A,B;C;x)  r (A) T (B)  ^ 


r(C)  r r(A+n)T(B+n)  x 
(A) T (B)  n£0  r(C+n)  n! 


(25) 


The  general  solution  of  equation  (23)  may  be  written,  for 
| x | < 1,  in  terms  of  the  original  variables  z and  w(z)  as 


- , , _ -zK/H„  , A n _ ft*  -z/H. 

w(z)  = D,e  F (A,B;C;  e ) + 


+ D2eZK/H  [^2  (KqI^2~  ) F(A-C+l,B-C+l;2-C;-p-^yi^y  e Z/H)  , 


(26) 


where  and  D2  are  arbitrary  constants. 


IV.  EIGENMODES  FOR  A RIGID  LOWER  BOUNDARY 


We  now  examine  modes  of  propagation  in  the  relatively  simple 


LThe  hypergeometric  nature  of  the  wave  equation  in  this  case  was 
noted  by  us  earlier  (Nye  and  Thomas  1974b) , and  also  independently 
by  Adam  (1974) . 
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case  of  an  isothermal  atmosphere  with  a uniform  horizontal 
magnetic  field,  bounded  from  below  by  a rigid  wall  at  z = 0. 

The  general  solution  (26)  is  subject  to  boundary  conditions 
at  z = 0 (x  = n2/B2 (Kz-n2 ) ) and  at  z = » (x  = 0) . 

As  a condition  at  z = °°,  we  require  that  the  total  energy 
of  the  perturbation  be  integrable  over  0 < z < The  magnetic 

energy  of  the  perturbation  is  proportional  to  the  square  of 


the  velocity.  Since 


1 im  F ( a , B ; Y ; x ) = 1 , 

x+0 


we  see  from  (26)  that  we  must  take  D2  = 0. 

The  second  boundary  condition  is  that  the  vertical  velocity 
vanish  at  the  rigid  wall,  i.e.,  w = 0 at  z = 0,  i.e.,  at 
x = n2/B2 (K2-ft2) . Provided 


B2 (K2-fi2) 


we  may  apply  this  condition  directly  to  (26)  with  D2  = 0 to 
obtain  the  dispersion  relation 


F ( A , B ; C ; g 3 (K-2_flg)  ) = 0 , (29] 

If,  however,  (28)  is  not  satisfied,  then  other  representations 
of  the  general  solution  (26),  valid  for  |x|  > 1 must  be  used  in 
order  to  apply  the  boundary  condition  at  z = 0.  For 


nf > i 

B2 (K2-n2)  - 
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the  interval  0 < x < nz/g2(Kz-ftz)  contains  the  regular  singular 
point  at  x = 1,  and  no  solution  that  satisfies  the  boundary 
conditions  and  is  also  regular  at  x = 1 is  to  be  expected. 

For  the  range 


ft2 

6 2 ( K 2 — 2 ) 


< - 1 , 


we  may  use  analytic  continuation  to  extend  the  general 
solution  (26).  The  analytic  continuation  of  (26)  with  D2  = 0, 
valid  for  x < - 1,  is  given  by 


wfz>  = D r~ZK/H  I T (C) T (B-A)  g2(fl2-K2)A  zA/H 

W(  ' Dl  | [ r (B) ] 2 l jp  J e 


r,«  t n.i  r,« . e2 (K2-n2)  z/h,  . r(or(A- 

F (A,  1-B,  1-B+A; ^ e ) + — f r ( mT 


x eZB/HF(B,l-A!l-A+B,iili^Sii  ez/H)j  . 


Thus,  for  the  range  of  parameters  (31) , the  dispersion  relation 


is  given  by 


r (B-A)  r32(02-K2)A  g2(K2-ft2)  , 

T r'(B ) F 1 SP ] f(a,i-b,i-b+a, ^ ) + 

r (a-b)  62(n2-K2),B  . r. 62 (K2-n2) . _ 

TtTaTTt  [ fil  1 ) - o . (33) 


The  dispersion  relation  (equation  (29)  or  (33))  has  been 
evaluated  for  various  values  of  the  nondimensional  parameter 
g2.  Examples  are  plotted  in  figures  1 and  2 (also  see  figure  4) 
The  curves  in  these  figures  represent  well-defined  eigenmodes 
of  wave  propagation  in  the  atmosphere.  These  curves  represent 
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trapped  waves  propagating  horizontally  in  the  waveguide  formed 
by  the  solid  boundary  below  and  the  exponentially  increasing 
Alfv£n  velocity  above. 

There  are  several  ways  of  interpreting  the  effect  of  , 

changing  g2  on  the  dispersion  relation.  First,  different 
values  of  g2  can  be  taken  to  represent  the  same  magnetic  field 
strength  and  the  same  density  at  z = 0,  but  different  atmospheric 
temperatures.  A second  interpretation  is  that  different  values 
of  g2  represent  the  same  magnetic  field  strength  and  the  same 
temperature,  but  different  densities.  This  is  equivalent  to 
placing  the  solid  lower  boundary  at  successively  higher  levels 
in  the  atmosphere  corresponding  to  larger  values  of  g2.  At 
each  higher  level  the  magnetic  restoring  force  becomes  more 
important  due  to  the  decreased  density,  while  the  compressible 
restoring  force  remains  the  same. 

The  third  interpretation  is  to  consider  changes  in  g2  to 
be  due  to  changes  in  the  magnetic  field  strength,  with  fixed 
values  of  temperature  and  density.  Since  the  sound  speed  and 
the  density  scale  height  then  do  not  change,  the  scales  for  the 
frequency  and  horizontal  wavenumber  are  the  same  in  each  case 
and  the  dispersion  diagrams  can  be  compared  directly.  From 
figures  1 and  2 it  can  be  seen  that  increasing  the  magnetic  field 
strength  (increasing  g2)  increases  the  cutoff  frequency.  As 
g2  increases,  the  slope  of  the  dispersion  curves,  and  hence  the 
group  velocity,  also  increases. 

Lowering  the  value  of  y to  represent  crudely  the  effect  of 


radiative  transfer  has  little  effect  on  the  nondimensional 


dispersion  diagrams.  However,  the  frequency  scaling  depends 
on  y as  w ~ (y)  2&,  while  the  wavenumber  scaling  is  independent 
of  y.  Therefore,  for  lower  y (lower  sound  speed)  the  wave 
oscillates  less  rapidly  and  the  phase  and  group  velocities  are 
correspondingly  lower. 

The  vertical  velocity  of  the  disturbance  (equation  (26)  or 
(32))  can  be  calculated  as  a function  of  height  for  any  point 
on  a dispersion  curve.  Figure  3 compares  the  lowest  mode  of 
oscillation  for  the  same  horizontal  wavenumber  but  different 
values  of  82  (i.e.,  different  magnetic  field  strengths),  and 

shows  that  for  increasing  magnetic  field  strength,  the  wave 
oscillates  more  rapidly  and  is  trapped  at  lower  levels  in  the 
atmosphere.  We  now  discuss  the  eigenmodes  given  by  the  dispersion 
relation  (29)  or  (33)  in  relation  to  an  observed  solar 
oscillation . 

V.  APPLICATION  TO  FLARE-INDUCED  CORONAL  WAVES 

On  September  20,  1963,  Moreton  and  Ramsey  (Moreton  1965) 
observed  a chromospheric  disturbance,  apparently  caused  by  the 
flash  phase  of  a flare,  propagate  at  a nearly  constant  velocity 
of  750  km  s ^ for  several  hundred  thousand  kilometers  across 
the  solar  disk.  Many  other  flare-induced  disturbances  have 
been  reported  (Moreton  1960,  Athay  and  Moreton  1961,  Dodson  and 
Hedeman  1968)  and  the  propagation  velocity  is  usually  on  the  order 


fc-’"  ■ - - ■■  ■ 
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of  1000  km  s ^ Dodson  and  Hedeman  (1968)  report  that  the 
width  of  the  disturbance  created  by  the  proton  flare  of 
August  28,  1966  was  greater  than  100,000  km. 


These  disturbances  could  not  have  been  propagating  solely 
in  the  chromosphere,  since  in  the  chromosphere  the  sound  speed 
is  only  of  the  order  of  20  km  s ^ and  the  Alfv4n  velocity  is 
only  of  the  order  of  50  km  s ^ Thus  a purely  chromospheric 
disturbance  would  have  created  a shock  wave  and  been  rapidly 
dissipated.  In  the  corona,  however,  both  the  Alfv^n  velocity  and 
the  sound  speed  are  an  order  of  magnitude  higher  than  in  the 
chromosphere  due  to  the  increased  temperature  and  decreased 
density.  It  has  been  proposed  (Meyer  1968;  Uchida  1968,  1970, 
1974;  Uchida  et  al.  1973)  that  the  disturbance  is  a magneto- 
hydrodynamic wave  propagating  in  the  low  corona  and  that  the 
motion  of  this  wave  at  the  corona-chromosphere  transition  region 
is  what  is  actually  observed.  There  is,  however,  no  general 
agreement  as  to  the  wavelength  of  the  disturbance  or  even  whether 
the  observed  disturbance  is  a single  wave  or  a wave  packet. 

Meyer  (1968)  studied  the  propagation  of  the  magneto-acoustic 
fast  mode  in  an  isothermal  corona  permeated  by  a uniform 
vertical  magnetic  field,  with  a rigid  lower  boundary  representing 
the  chromosphere-corona  transition  region.  He  found  eigenmodes 
with  nearly  constant  horizontal  group  velocity.  Equating  the 
group  velocity  to  the  observed  propagation  velocity,  Meyer 
found  that  for  a horizontal  wavelength  on  the  order  of  100,000  km, 
the  coronal  magnetic  field  must  be  approximately  6 gauss,  a 
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n — pm 


1 


reasonable  average  value. 

Uchida,  in  a series  of  papers  (Uchida  1968,  1970,  1974; 

Uchida  et.  al.  1973),  studied  the  propagation  of  short  wave- 
length (-5000  km)  disturbances  in  various  realistic  coronal 
models.  Using  a ray-tracing  technique,  he  obtained  horizontal 
and  vertical  refractions  in  close  agreement  with  the  observed 
waves . 

Although  the  magnetic  field  structure  of  the  corona  is 
quite  complicated,  the  field  changes  fairly  slowly  and  there 
are  probably  regions  of  nearly  uniform  field  with  almost  any 
orientation.  We  study  the  case  of  a uniform  horizontal  field 
in  connection  with  the  coronal  wave  problem  only  as  a means  of 
understanding  the  mechanism  of  wave  propagation  for  waves  of 
arbitrary  wavelength.  Our  model  supplements  Meyer's  (1968)  work 
by  considering  the  case  of  a uniform  horizontal  magnetic  field, 
and  by  including  effects  of  gravity  and  stratification.  As  in 
Meyer's  model,  we  use  the  rigid  lower  boundary  to  represent 
upward  reflection  from  the  chromosphere-corona  transition  layer. 

We  have  evaluated  our  solution  for  a temperature  of 
1.6  x 106  °K  and  g2  = 10,  which  is  fairly  typical  of  the  base 
of  the  corona.  These  values  correspond  to  a sound  speed  of 
180  km  s ^ and  a density  scale  height  of  approximately  71,000  km. 
Figure  4 shows  the  dimensional  dispersion  relation  for  these 
parameters.  We  see  that  the  dispersion  curves  are  nearly  straight, 
which  means  that  these  modes  have  very  little  dispersion  and  will 
propagate  for  great  distances  with  little  change  in  character. 
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The  first  three  modes  have  been  calculated  for  a wavelength 
of  75,000  km  in  figure  5.  The  first  mode  has  nearly  zero 
vertical  velocity  above  two  scale  heights  and  is  therefore 
trapped  in  the  low  corona. 

The  phase  velocity  and  group  velocity  of  the  first  mode 
have  been  plotted  as  a function  of  horizontal  wavenumber  in 
figure  6.  For  any  wavelength  of  100,000  km  or  less,  the 
group  velocity  is  nearly  constant  at  about  610  km  s \ Since 
the  energy  of  a disturbance  propagates  at  the  group  velocity, 
it  is  not  important  which  specific  wavelength,  or  spectrum 
of  wavelengths,  receives  energy  from  the  flare.  The  energy 
at  all  wavelengths  will  propagate  together  as  a wave  packet 
near  the  lower  coronal  boundary. 

The  present  model  is  not  proposed  as  a realistic  model  of 
the  solar  corona,  although  it  may  be  fairly  accurate  over 
certain  regions.  No  attempt  has  been  made  to  include  the 
effects  of  horizontal  variations.  The  value  of  the  model  is 
that  a mechanism  for  wave  propagation  can  be  studied  for 
arbitrary  wavelengths.  These  results  close  the  gap  between 
the  short  wavelength  ray-tracing  theory  and  the  long  wavelength, 
vertical  field  case.  We  show  that  the  question  of  wavelength 
is  not  particularly  important  since  the  group  velocity  of  the 
trapped  modes  is  essentially  independent  of  wavelength.  For  the 
relatively  large  value  of  B2  (=10),  the  wave  modes  are  basically 


the  magneto-acoustic  fast  modes  (studied  by  Meyer  and  by  Uchida) 
modified  by  gravity.  For  an  inclined  magnetic  field,  there  must 
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also  be  trapped  modes  of  propagation  involving  a coupling  of 
the  present  modes  and  the  type  of  mode  studied  by  Meyer  for  a 
vertical  field.  The  present  results,  taken  with  those  of 
Meyer  and  the  work  of  Uchida,  present  a consistent  picture  of 
flare-induced  coronal  waves  as  guided  magneto-atmospheric 
waves . 
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FIGURE  CAPTIONS 


Fig.  1 


Fig.  2 
Fig . 3 


Fig . 4 


Fig.  5 

Fig.  6 


Dispersion  diagrams  (nondimensional  frequency  versus 
nondimensional  horizontal  wavenumber)  for  an  isothermal 
atmosphere  with  a uniform  horizontal  magnetic  field 
and  a solid  lower  boundary,  with  y = 5/3,  B2  = 10“ 5 , 
and  B2  = 10“ 4 . The  curves  represent  eigenmodes  and  the 
crosses  indicate  the  computed  points.  This  figure 
should  not  be  confused  with  a diagnostic  diagram  for  an 
atmosphere  with  constant  parameters  (e.g.,  Yu  1965). 

Same  as  Fig.  1,  but  with  B2  = 0.5. 

Direct  comparison  of  the  first  mode  of  oscillation 
of  the  atmosphere  for  the  same  nondimensional  horizontal 
wavenumber  (K  = 0.25)  but  for  two  values  of  B2 
(B2  = 10'3,  B2  = 10“ 1 ) . 

Dimensional  dispersion  diagram  for  B2  - 10,  y = 5/3, 
and  T = 1.6  x 106oK.  The  curves  represent  eigenmodes 
of  the  corona  which  are  trapped  by  the  increasing 
Alfvdn  velocity  with  height. 

First  three  modes  of  oscillation  of  the  model  corona 
for  B2  = 10,  y = 5/3.  The  vertical  velocities  have 
been  normalized  to  maximum  value  unity. 

Phase  velocity  Vp  and  group  velocity  v^  of  the  first 
mode  of  coronal  oscillation  plotted  as  a function  of 
horizontal  wavenumber. 
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ABSTRACT 


F 


A simple  two-layer  model  of  a sunspot  penumbra  is  used 
to  study  the  mode  of  running  penumbral  waves.  Exact  solutions 
of  the  .linearized  wave  equation,  not  limited  to  the  small- 
wavelength  approximation,  are  employed  in  each  layer.  The 
lowest  "plus"  eigenmode  of  magneto-atmospheric  waves  in  the 
model  penumbra  is  in  good  agreement  with  observations  of 
running  penumbral  waves.  The  results  indicate  that  running 
penumbral  waves  should  be  observable  in  a photospheric  spectral 
line. 
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I. 


INTRODUCTION 


In  Paper  I of  this  series  (Nye  and  Thomas  1975)  we 
presented  an  exact  analytical  solution  for  magneto-atmospheric 
waves  in  the  case  of  an  isothermal  atmosphere  with  a uniform 
horizontal  magnetic  field.  In  the  present  paper  we  apply  this 
solution  to  a simple  two-layer  model  penumbra  in  order  to  study 
the  mode  of  running  penumbral  waves. 

Running  penumbral  waves  (Zirin  and  Stein  1972;  Giovanelli 
1972,  1974;  Moore  and  Tang  1975)  are  good  examples  of  magneto- 
atmospheric  waves.  These  waves  propagate  radially  outward 
across  sunspot  penumbrae,  with  predominantly  vertical  motions 
in  Ha.  The  observed  range  of  frequency  and  propagation  speed 
is  fairly  well  established  (see  discussion  in  section  IV) . 

Moore  (1973)  has  concluded  that  the  source  of  excitation  of 
the  penumbral  waves  is  overstable  convection  in  the  low  umbra. 

In  an  earlier  paper  (Nye  and  Thomas  1974  (NT) ) we  studied  the 
mode  of  propagation  of  penumbral  waves  on  the  basis  of  a piece- 
wise  linear  model  of  the  vertical  structure  of  a typical  sunspot 
penumbra.  We  found  the  penumbral  waves  to  be  magneto-atmospheric 
waves  (of  the  "plus"  type)  that  are  vertically  trapped  at  photo- 
spheric  levels.  This  trapping  is  primarily  due  to  the  increasing 
sound  speed  with  depth  into  the  convection  zone  and  the  increasing 
Alfv£n  velocity  with  height  into  the  chromosphere. 

Here  we  extend  our  earlier  work  by  computing  actual  eigen- 
modes  of  propagation  for  a somewhat  simpler  model  penumbra,  which 
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nevertheless  retains  the  essential  features.  The  properties  of 
the  lowest  mode  of  propagation  of  the  model  penumbra  turn  out 
to  be  in  good  agreement  with  observations,  and  give  some  useful 
clues  for  further  observation  of  running  penumbral  waves. 

II.  THE  TWO-LAYER  PENUMBRAL  MODEL 

The  entire  penumbral  model  consists  of  a compressible, 
inviscid,  perfectly-conducting,  stratified  perfect  gas  subject 
to  a constant  acceleration  of  gravity  g (=  0.274  km  s-2)  in  the 
negative  z-direction.  The  upper  layer  is  isothermal  and  is 
permeated  by  a uniform  horizontal  magnetic  field,  which  yields 
an  Alfven  velocity  that  increases  exponentially  with  height  due 
to  the  decreasing  density.  An  exact  solution  of  the  linearized 
propagation  equation  for  this  case  was  given  in  Paper  I. 

This  upper  layer  is  a suitable  model  of  the  penumbral  photo- 
sphere and  chromosphere,  where  observed  penumbral  magnetic  fields 
are  very  nearly  horizontal  (Nishi  and  Makita  1973)  and  decrease 
slowly  with  height  (Bray  and  Loughhead  1964).  The  scale  height 
for  variation  of  the  magnetic  field  is  very  large  compared  to 
the  density  scale  height,  so  the  assumption  of  a uniform  horizontal 
magnetic  field  is  reasonable.  Our  earlier  calculations  (NT)  showed 
that  running  penumbral  waves  are  trapped  at  photospheric  levels, 
so  that  the  increasing  sound  speed  in  the  upper  chromosphere  has 
little  effect  on  the  trapping.  Taking  the  upper  layer  to  be  iso- 
thermal is  therefore  also  a reasonable  assumption. 


The  vertical  distributions  of  the  sound  speed  and  Alfv^n 
velocity  for  the  two-layer  model  penumbra  are  shown  in  Figure  1. 
Subscripts  1 and  2 denote  quantities  in  the  upper  and  lower 
layers,  respectively,  and  the  subscript  0 refers  to  quantities 
evaluated  at  z = 0.  The  sound  speed,  density  scale  height,  and 
Alfven  velocity  in  the  upper  layer  are  given  by 


and 


Cj  = YjRTj  = const.  , 
2 

C1 

H,  = — — = const.  , 

1 Y,g 


2 , , 2 Z/Hi 

VA(Z)  = V e 


(1) 

(2) 

(3) 


The  lower  layer  of  the  penumbral  model  (layer  2,  figure  1) 

is  adiabatic  with  no  magnetic  field.  The  temperature  decreases 

with  height  (increases  with  depth)  at  the  adiabatic  lapse  rate, 

(dT/dz)  = -g/c  , and  thus  this  layer  is  neutrally  stable.  The 
s p 

actual  temperature  distribution  in  the  convection  zone  below  a 
penumbra  is  probably  very  nearly  adiabatic,  except  for  a thin 
superadiabatic  layer  just  beneath  the  photosphere  that  we 
neglect  here.  There  is  no  magnetic  field  in  this  layer  since 
we  assume  that  the  penumbral  magnetic  field  lies  over  the  con- 
vection zone.  The  sound  speed  squared  and  local  density  scale 
height  each  increase  linearly  with  depth  in  the  lower  layer, 
their  functional  forms  being 

c2  (z)  = c2  0 - g (y2-1)  z (4) 

and 
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(5) 


H 2 (z)  = H2  0 - ( Y 2 — 1 ) z 


The  corresponding  density  distribution  is 

P2(z)  = P20  [1  - (Y2-Dz/H20]1/(Y2-^  . (6) 


At  the  interface  between  the  two  layers  (z=0) , we  require 
the  undisturbed  density  to  be  continuous  to  avoid  introducing 
interfacial  gravity  waves  and  wave,  reflections;  therefore, 

P10  = p . In  the  unperturbed  penumbra,  there  must  be  pressure 
equilibrium  at  the  interface;  that  is,  the  gas  pressure  in  the 
lower  layer  at  z=0  must  equal  the  sum  of  the  gas  pressure  and 
the  magnetic  pressure  in  the  upper  layer  at  z=0.  Therefore, 
the  gas  pressure  is  greater  in  layer  2 than  in  layer  1,  and 
since  the  density  is  continuous  across  the  interface,  the  tem- 
perature is  greater  in  layer  2.  This  may  be  expressed  in  terms 
of  the  sound  speeds  and  the  Alfv^n  velocity  at  z=0  as 


c 


2 

2 0 


(7) 


We  now  turn  to  the  problem  of  computing  eigenmodes  of 
magneto-atmospheric  waves  in  this  penumbral  model. 


III.  ANALYSIS 

Consider  first  the  behavior  of  small  adiabatic  perturbations 
in  the  lower  layer  (layer  2) , an  adiabatic  atmosphere  without 
magnetic  field.  Leibacher  (1971)  solved  this  problem  in  his  study 
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of  oscillations  of  the  quiet  photosphere.  Here  we  take  a 
slightly  different  approach  than  his,  using  other  transformations 
which  yield  a different  form  of  the  propagation  equation.  For 
vanishing  magnetic  field  (B^  = B = 0) , the  vector  equation  for 
the  perturbation  velocity  (equation  (7) , Paper  I)  becomes 

3 2u 

9 jy  = c^v<(>  + (y  2-l ) <}>g  + V(u2*g_)  , (8) 

where  <J>  = V*u2- 

We  assume  that  the  perturbation  velocity  has  the  form 
u 2 = u2(z)e  x , with  propagation  in  the  x-direction  (k  = 0) . 

This  implies  that  v2  = 0 and  u2  = u2(z)  = (u2 (z) ,0,w2 (z) ) . 

All  other  perturbation  quantities  are  represented  in  a similar 
manner,  with  a caret  denoting  the  z-dependent  amplitude  in 
each  case.  From  the  two  components  of  equation  (8)  and  the 
definition  of  <j> , we  obtain  the  following  relation: 


w2  (z) 


g [Y2~c2kx/ u)2]  <J>  - c2 
[w2-g 2 k£/ cjo  2 ] 


d<j> 

dz 


(9) 


The  pressure  perturbation  can  also  be  written  in  terms  of  <f> , 
using  the  continuity  and  energy  relations  (see  Paper  I) : 


P2  = 


iP2C2  j 

“ t 


2 ( [w2-g  2y2/ c * ] $ + 


[ w 2 -g 2 k^/ oj  2 ] 


* all 


(10) 


Upon  substitution  of  equation  (9)  into  the  z-component  of 
equation  (8),  we  obtain  the  following  second-order  differential 
equation  for  4> : 
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M ( a , b ; Y ) 


r (b)  " r (a+n)  Yn 
r (a)  n£0  r(b+n)  n!  * 


The  general  solution  of  equation  (14)  is 


(16) 


(17) 


where  D3  and  D 4 are  arbitrary  constants  and  U(a,b;Y)  can  be 
written  in  terms  of  Rummer's  functions.  To  insure  finite  total 
perturbation  energy,  we  must  require  that  the  vertical  velocity 
of  the  perturbation  vanish  as  z -*■  -<*>  or  as  Y -*  +«>.  This  in 
turn  requires  that  D3  = 0.  The  solution  for  J as  a function 
of  z in  the  lower  layer  is  then 


The  vertical  velocity  and  pressure  perturbation  are  given  in 
terms  of  $ by  equations  (9)  and  (10). 

The  form  of  the  solution  in  the  upper  layer  has  been  given 
in  Paper  I,  and  we  shall  not  repeat  the  analysis  here.  The 
general  form  for  the  vertical  velocity  in  the  upper  layer  that 
gives  finite  total  perturbation  energy  is  given  in  terms  of 
hypergeometric  functions  by  either 


w,(z)  = 


F ( A, B ; C ; 


8 2 0<2-n2) 


,-z/Hi 


for 


or  else  by 


X 


Wj  (2) 


D 1 e 


-zK,/H 


( r (c)  r (b-a) 

reM^-K2)] 

[ r (B)  ] 2 

L 

eZA/H, 


F ( A , 1-B ; 1-B+A ; 


r (o  r (a-b)  re2  f-K i ) 

ir (A) ] 2 ft2 


B 


zB/H  . 
e 1 


X 


F (B , 1-A; 1-A+B ; 


B2 ( K2  -ft  2 ) 


ft? 


(21) 


for 


Y 

B2  (K 2 -ft 2 ) > 

Here 

A + B = C = 2K 

AB  = fij  + K,  + 
B2  = v 2 /c 2 , 


ez/Hi  . (22) 


+ 1 , 

(23) 

. 2 

IilL 
y i 'Qi  ' 

(24) 

(25) 

and  ft j and  K:  are  nondimensional  frequency  and  wavenumber  defined 
as  in  (12),  except  scaling  with  Cj  and  Hj. 

The  pressure  perturbation  in  the  upper  layer  consists  of  the 
sum  of  the  gas  pressure  perturbation  and  the  magnetic  pressure 
perturbation.  The  gas  pressure  perturbation  can  be  expressed 
in  terms  of  the  vertical  velocity  and  its  derivative  as 


Pi  = 


-ip  ,C  1 CO 


- Ar  W . 
2 1 dz  c l_  1 


(oi2-c^kM 


(26) 


The  magnetic  pressure  perturbation  pm  is  found  from  the 
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ifill 





linearization  of 


p + p = 
'm 


(B 0+B)  • (B 0+B) 

8 7T 


where  pm  is  the  unperturbed  magnetic  pressure.  The  components 
of  the  perturbed  magnetic  field  are  determined  by  the  linearized 
induction  equation  (equation  (4)  of  Paper  I).  We  find 


~ . VA^ 1 dw, 

^m  1 to  d z 


(27) 


The  total  perturbed  pressure  p in  layer  1 is  expressed  in 

1 i 

terms  of  the  density  and  the  Alfven  velocity  as 


i^{  I (cj+vi)»*-c;vjk;]^- w2gwl}  / (u>’-c?k£)  . 


(28; 


IV.  EIGEN MODES  AND  RUNNING  PENUMBRAL  WAVES 

We  now  have  expressions  for  the  vertical  velocity  and  the 
pressure  perturbation  in  each  layer  of  the  penumbral  model,  such 
that  the  total  oerturbation  energy  is  finite.  The  remaining 
conditions  are  the  matching  of  the  vertical  velocity  and  the 
perturbed  pressure  at  the  interface  z = 0. 

The  scaling  of  frequency  and  wavenumber  was  done  separately 
for  each  layer  in  order  to  simplify  the  propagation  equation  as 
much  as  possible  in  each  case.  In  matching  across  the  interface, 
we  need  the  following  relations  between  parameters  in  the  two 
layers : 

»20  = Yad  + , (29) 
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(30) 


g 2 y 

Ql  = YjY2(1  + 1-rL)^  ' 

K 2 = y*  (1  + — 2-M  K2  , (31) 

and  equation  (7) . 

Continuity  of  the  vertical  velocity  across  z = 0 requires, 
after  normalization,  that 

Wj  (0)  = w2  (0)'  = 1 . (32) 

This  condition  fixes  the  values  of  the  coefficients  Dj  in  equation 
(19)  or  (21)  and  04  in  equation  (18).  The  remaining  condition, 
the  continuity  of  the  perturbed  pressure,  requires  that  we 
equate  (10)  and  (28).  This  then  leads  to  the  nondimensional 
condition 

I t(i+e*)B;-s*K;i*a|rSf  H,(0)  j l (ii,T,«(0)  + ^ ) 

I (TTpiTf)  j y,  j TQpK^7nf)  j ' {33) 

where  0 = $H20.  Equation  (33)  is  only  satisfied  by  particular 
values  of  frequency  and  wavenumber,  and  gives  the  dispersion 
relation  for  eigenmodes  of  oscillation  in  the  penumbral  model. 

In  order  to  evaluate  the  dispersion  relation  (33)  , the  free 
parameters  62,  y1 , and  y2  must  be  specified,  which  then  effectively 
determines  the  properties  of  the  model.  Although  it  is  possible 
to  evaluate  (33)  for  different  values  of  y in  each  layer  (for 
example,  a lower  value  of  y could  be  taken  to  represent  radiative 
transfer  in  the  upper  layer) , we  chose  the  usual  value  of  5/3  for 
both  layers.  The  dispersion  relation  (33)  was  solved  numerically 
by  inserting  values  of  and  then  computing  and  comparing  the  two 
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sides  of  (33)  for  small  increments  in  ft  . 

Figure  2 shows  the  first  several  eigenmodes  of  the  two- 
layer  penumbral  model  for  B2  = 0.5  and  Yj  = Y2  = 5/3.  The 
value  of  B2  was  chosen  to  represent  a typical  penumbra  and  is 
slightly  less  than  the  value  of  S!  at  z = 0 in  our  earlier 
penumbral  model  (NT) . Here  we  have  classified  the  eigenmodes 
as  "plus"  or  "minus"  modes,  following  the  terminology  used  in 
the  case  of  an  atmosphere  with  constant  sound  speed,  density 
scale  height,  and  Alfven  velocity  (see  KcLellan  and  Winterberg 
1968  and  NT).  The  plus  modes  all  lie  above  the  upper  dashed  line 
ft  = K in  figure  2,  which  corresponds  to  u = ci^’x  Lam^  mode). 

The  minus  mode  (there  is  only  one  in  this  case)  lies  below  the 
lower  dashed  line  ft  = BK  which  corresponds  to  w = vQkx.  There 
are  no  eigenmodes  in  the  region  between  the  dashed  lines.  This 
classification  of  plus  or  minus  modes  refers  here  to  the  character 
of  the  eigenmode  in  the  upper  layer;  in  the  lower  layer,  all  of 
the  modes  have  the  character  of  acoustic  waves  (no  magnetic 
field  and  no  buoyancy) . 

Observational  data  on  penumbral  waves  are  included  in  figure 
2 for  comparison.  The  most  commonly  reported  observational 
quantities  are  the  period  and  the  horizontal  phase  velocity, 
although  they  are  not  always  measured  simultaneously.  Giovanelli 
(1974)  reports  a typical  phase  velocity  of  running  penumbral  waves 
of  15  km  s'1  and  typical  periods  in  the  range  of  180-240  s.  He 
did  report  phase  velocities  of  up  to  21  km  s-1,  however,  and  gave 
specific  periods  and  wavelengths  for  four  sunspots  (denoted  by 
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crosses  in  figure  2) . The  data  of  Beckers  and  Schultz  (1972) 
appear  to  indicate  a penumbral  oscillation  period  of  255  s. 

Moore  and  Tang  (1975)  observed  penunbral  waves  with  period 
270±10  s in  a single  sunspot.  Zirin  and  Stein  (1972)  state 
that  the  periods  of  penumbral  waves  in  twenty  sunspots  were 
almost  all  between  240  s and  300  s,  and  the  measured  horizontal 
phase  velocity  of  9.4  km  s-1  in  one  spot  was  more  or  less  the 
same  in  other  spots  even  when  the  period  varied. 

The  quadrangle  in  figure  2 represents  the  range  of  observa- 
tions: periods  from  180  s to  300  s,  and  phase  velocities  from 

9.4  km  s-1  to  21  km  s-1,  with  a dashed  line  at  15  km  s-1  to 
indicate  the  value  that  Giovanelli  considers  typical.  The 
first  plus  mode  of  the  penumbral  model  passes  through  this 
quadrangle.  Although  the  particular  eigenmode  of  oscillation  of 
the  penumbra  is  determined  by  the  excitation,  and  little  is  known 
about  the  excitation,  the  present  results  indicate  that  it  is  the 
first  plus  eigenmode  that  is  being  excited.  This  agrees  with  our 
earlier  conclusions  (NT) . 

The  value  f32  = 0.5  used  in  figure  2 was  chosen  to  represent 
a typical  penumbra.  In  figure  3 the  effect  of  changes  in  (32  on 
the  first  plus  mode  is  shown  for  a range  of  82  of  two  orders  of 
magnitude  (0.05,  0.5,  5.0).  This  constitutes  a reasonable  set 
of  limits  on  8 2 for  penumbral  conditions,  and  is  obtained  by 
looking  at  the  normal  variation  of  B (factor  of  4) , p (factor  of 
4),  and  c2  (factor  of  1.5)  expected  in  different  penumbrae. 

We  see  that  for  any  reasonable  value  of  82,  the  penumbral  model 
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has  a first  plus  mode  within  the  range  of  observations. 


The  vertical  distributions  of  velocity  and  kinetic  energy 
of  the  first  plus  mode  (for  B2  = 0.5)  are  shown  in  figure  4 
for  Kj  = 0.1995  and  fij  = 0.3628,  corresponding  to  a horizontal 
wavelength  \ = 3,000  km  and  period  t = 250  s for  c2  = 43.5  km2s“2. 
Here  the  nondimensional  height  is  scaled  everywhere  by  the 
density  scale  height  in  the  upper  layer,  Hj.  The  velocity 
distribution  is  fairly  symmetric  with  the  maximum  amplitude 
occurring  slightly  above  z = 0 in  the  penumbral  photosphere. 

The  kinetic  energy,  on  the  other  hand,  is  almost  entirely  trapped 
in  the  lower  layer  (convection  zone)  with  maximum  energy  just 
below  the  interface.  The  velocity  amplitude  decays  slowly  with 
height  with  a value  of  more  than  25%  of  the  maximum  amplitude 
at  a distance  of  eight  scale  heights  above  the  level  of  that 
maximum. 

There  is  a discrepancy  between  the  height  of  maximum 
velocity  predicted  here  (z  ~ 100  km)  and  that  predicted  in  our 
earlier  paper  (NT,  z ~ 1000  km) . Here  the  Alfv£n  velocity 
increases  exponentially  with  height  above  the  photosphere, 
whereas  in  NT  it  increased  linearly.  Thus,  the  downward  refrac- 
tion of  waves  is  much  stronger  in  the  present  model.  The  actual 
situation  is  probably  somewhere  between  these  two  cases.  In 
either  case,  the  wave  energy  lies  mostly  below  the  height  of 
maximum  vertical  velocity,  in  the  convection  zone  and  low  photo- 
sphere. 
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V.  CONCLUSIONS 


The  present  results,  taken  together  with  our  earlier  work 
(NT),  indicate  that  running  penumbral  waves  should  be  identified 
with  the  lowest  plus  mode  of  trapped  magneto-atmospheric  waves 
in  the  penumbra.  The  vertical  trapping  is  primarily  due  to  the 
increasing  Alfven  velocity  up  into  the  chromosphere  and  the  in- 
creasing sound  speed  down  into  the  convection  zone.  Most  of 
the  energy  of  the  penumbral  waves  lies  in  the  convection  zone 
and  low  photosphere,  at  the  same  level  as  the  expected  source 
of  excitation  (umbral  oscillatory  convection) . The  maximum  wave 
amplitude  occurs  somewhat  higher. 

The  results  also  indicate  that  penumbral  waves  should  be 
observable  in  a photospheric  spectral  line  (see  figure  4)  as 
well  as  in  Ha.  There  is  some  indication  of  this  in  the  obser- 
vation of  Beckers  and  Schultz  (1972) . Their  data  show  a 255  s 
period  oscillation  in  the  penumbra  of  one  sunspot  observed  in 
a photospheric  line.  They  present  contours  of  vertical  velocity 
as  a function  of  horizontal  position  and  time  (their  figure  1) 
in  which  one  may  note  horizontal  propagation  outw’ard  across 
the  penumbra  at  about  the  right  phase  speed.  We  plan  further 
observations  in  a search  for  penumbral  waves  in  the  photosphere. 
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gradient  and  no  magnetic  field. 

Fig.  2 The  first  five  plus  eigenmodes  and  the  only  minus 

eigenmode  of  the  two-layer  penumbral  model  for  B2  = 0.5, 
Yj  = Y2  = 5/3.  The  quadrangle  represents  the  range  of 
observational  data  (see  text).  The  crosses  correspond 
to  particular  observations  (Giovanelli  1974). 

Fig.  3 The  first  plus  eigenmode  evaluated  for  extreme  values 
of  B2  for  sunspot  penumbrae  (see  text):  82rax  = 5.0, 
B2mjn  = 0.05.  The  quadrangle  represents  the  range  of 
observational  data.  The  crosses  correspond  to 
particular  observations  (Giovanelli  1974). 

Fig.  4 The  distribution  of  vertical  velocity  and  kinetic  energy 
of  the  first  plus  eigenmode  with  nondimensional  height 
z/K j . The  velocity  and  energy  are  each  normalized  to 
value  unity  at  the  interface  z/H,  = 0. 


